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Abstract 


State Departments of Transportation (DOTs) need road ice condition forecasts in order to 
arrange anti-icing/de-icing activities and to enhance public safety and awareness. Icy road 
fatalities account for 3.6 times more deaths than all other weather hazards combined (USDOT, 
2001). The relative risk of vehicle crashes significantly increases given the presence of winter 
precipitation (Black & Mote 2015; Yu et al. 2019). Applying anti-icing materials on the roads, 
has been reported as having a cost reduction of 30%-90% compared to de-icing after road ice 
forms (Du et al. 2019). In order to help improve public safety this research was divided into two 
focus areas, Task 1 Development of IcyRoad Model and Task 2 Validation of the Ice Formula. 

Task 1 Development of IcyRoad Model - For Montana DOT (MDT), a statewide IcyRoad 
model has been developed, tested, and refined during winter 2020-2021 using Montana Road 
Weather Information System (RWIS) measurements (Ewan & Al-Kaisy, 2017). This product is 
demonstrated at https://sg-weather.com/INTERNAL/Icyroad3/. The MDT IcyRoad project 
identified the most accurate forecast algorithm for Montana (IcyRoad3), which assimilated 
RWIS hourly historical observations into a weather model forecasting temperature, clouds, 
rainfall information, and the SpringGem road-surface physical scheme for road surface layer 
turbulence and heat transportation. IcyRoad project achieved several milestones, including 
developing an online automatic RWIS validation tool, and developing an in-depth understanding 
of the formation of road ice. The overall hourly road ice forecast accuracy improved from an 
averaged 62% to 82%. 

Though progress on road ice formation mechanisms and forecasting capabilities has been 
made, research on black ice was very limited due to a lack of black ice identification and funding 
availability. More development on model forecast schemes for black ice formation is also 
needed. 

Task 2 Validation of the Ice Formula — Validating and refining the IcyRoad scientific 
algorithm, in particular the black ice algorithm was conducted during the project period spanning 
July 2020 - July 2021, progress was made investigating technologies, particularly the use of 
hyperspectral imagery, to identify black ice formation on road surfaces. Black ice mechanisms 
are one of the most meteorologically challenging processes to model and more observations are 
needed to understand black ice forming meteorological conditions, locations, and timing of 
events in order to forecast it well. 

An ice sheet, thin or thick, has spectral reflectance (rd) significantly different at near 
infrared bands. By using the combination of visible bands and near-infrared bands (NIR), a thin 
ice-covered surface may be differentiated from road surfaces. Theoretically, for snow and ice, ra 
reduces from RED to NIR and thus @ is negative; while for asphalt surfaces, rA increases from 
RED to NIR and thus has a positive 8. As a result, icy roads can theoretically be detected using a 
combination of Red and NIR bands. (Road Ice Index 0 = (NIR-RED)/(NIR+RED)). 

There is evidence of hyperspectral cameras with visible and near infrared capabilities (0.4 
— 1 um) being used to identify road defects (Abdellatif et al. 2019) but no one has yet attempted 
to distinguish road surface conditions using hyperspectral cameras. Therefore, the team has been 
developing a drone-based remote sensing technology to detect black ice via hyperspectral camera 
launched on an unattended aerial vehicle (UAV). This project used several hyperspectral 
cameras, a Pika L that covers the spectral range from 0.4 — 1 um and a Pika NIR-320 that covers 
the spectral range from .9-1.7 um. Additionally, a spectrometer with a range of .35—2.5 um was 
also utilized. 


Results within a controlled laboratory setting indicate promise for detecting ice versus 
dry asphalt using hyperspectral cameras but discerning black ice from other ice is not promising 
using the ice index formula with this current dataset. Another analysis method using a common 
method for processing and interpreting hyperspectral data with principal component analysis 
(PCA) (Demar et al., 2012) was also explored. This study demonstrated that the PCA technique 
is capable of identifying, water, and snow in comparison to dry asphalt. 

Stepping up to a UAS based collection platform there is a noticeable decrease in data 
resolution, due to sensor to target distance, subsequent accuracy decreases compared to the 
laboratory and field-based data. 


Task 1 Development of IcyRoad Model 
SECTION I 


Executive Summary 


SpringGem Weather Information, LLC is the sub-awardee of “Validation and 
Refinement Of Icy Road Forecast and Alert (IcyRoad) System Using MDT RWIS Sites” project 
sponsored by the Montana Department of Transportation (MDT) to the University of Montana 
(UM) (Principle Investigator (PI), Jennifer Fowler). During the project period spanning July 15, 
2020 - July 15, 2021, encouraging progress has been made by SpringGem Weather Information 
and UM in data collection, processing, and analysis. This task report outlines the milestones, 
scientific understandings, progress achieved, current status of IcyRoad algorithm development, 
and linkage with the National Science Foundation (NSF) and other Department of Transportation 
(DOT) IcyRoad projects. 

The first generation of the Icy Road Detection and Alert (1.e., IcyRoad) algorithm, funded 
by NSF SBIR Phase 1, is shown on a web application at https://sg-weather.com/weather/ as well 
as in a mobile application named “IcyRoad Alert”, together with value-added highway camera 
information http://www.sgweatherinfo.com/traffic_and a digital weather forecast. Based on 
Internet-of-Things technologies, IcyRoad integrates the numerical weather forecasts, road 
physical scheme developed in SpringGem Weather Information (Gohil & Jin 2019; Zhang et al. 
2019), ground observations, commercial cloud computing, data visualization and mining, and 
remote sensing technology to predict the road ice conditions for any road across the US with a 
24-hour lead time. The next generation of IcyRoad will offer a lead time of 72 hours to 8-day 
forecast and with improved forecast accuracy. 

The UM is working with SpringGem Weather Information through this project to validate 
and refine the IcyRoad scientific algorithm, in particular the black ice algorithm, using Montana 
Department of Transportation (MDT) Road Weather Information System (RWIS) data and UM’s 
drone-based ice detection technology particularly the use of hyperspectral imagery, to identify 
black ice formation on road surfaces. 


Introduction 


State Departments of Transportation (DOTs) need road ice condition forecasts in order to 
arrange anti-icing/de-icing activities and to enhance public awareness and safety. Applying anti- 
icing materials on the roads, a proactive process which involves applying anti-icing materials on 
the road surface before ice forms, has been reported as having a cost reduction of 30%-90% 
compared to de-icing after road ice forms (Du et al. 2019). For example, the “Winter Parking Lot 
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and Sidewalk Maintenance Manual” used in Minnesota states that anti-icing treatments applied 
in advance of stormy weather often cost an order of magnitude less than de-icing operations 
implemented after a storm (www.sima.org). Therefore, knowing when and where icy road 
conditions will occur is highly desired by state DOTs for both road safety and winter road 
maintenance material cost reduction. 

For Montana DOT (MDT), a statewide IcyRoad model has been developed, tested, and 
refined during winter 2020-2021 using Montana Road Weather Information System (RWIS) 
measurements (Ewan & Al-Kaisy, 2017). This product is demonstrated at https://sg- 
weather.com/INTERNAL/Icyroad3/!. The MDT IcyRoad project is an additional research effort 
to an NSF/SBIR-supported IcyRoad Development project by validating a set of IcyRoad 
forecasting algorithms using RWIS data. The MDT IcyRoad project identified the most accurate 
IcyRoad forecast algorithm for Montana (i.e., IcyRoad3)', which assimilated RWIS hourly 
historical observations into a weather model forecasting temperature, clouds, rainfall 
information, and the SpringGem road-surface physical scheme for road surface layer turbulence 
and heat transportation. 

Beyond weather conditions, road texture and geographic conditions affect road ice 
formation. For example, when temperature and water conditions are optimal, bridges or road 
shoulders (with even a slight slope) can form ice much more easily due to the micro-scale 
turbulences near the ground. Additionally, urban development also affects local convection, 
clouds condensation nuclei, and snowfall and thus road ice (Jin & Shepherd 2005; Jin et al. 2005; 
Jin & Shepherd 2008; Mote 2008; Jin & Dickinson 2010). Ground measurements from RWIS 
(Ewan & Al-Kaisy 2017; Johnson & Shepherd 2018) have been analyzed to better understand ice 
and black ice formation mechanisms across US roadways. 

Icy road fatalities account for 3.6 times more deaths than all other weather hazards 
combined (USDOT, 2001). In the decade from 2002 to 2012, winter weather conditions 
contributed to 540,000 car collisions (10% of collisions annually), 150,000 car-related injuries, 
and 1,900 car-related fatalities, according to the US Federal Highway Administration (Harris 
2018). The relative risk of vehicle crashes significantly increases given the presence of winter 
precipitation (Black & Mote 2015; Yu et al. 2019). 


Progress, Limits, and Results 


The MDT IcyRoad project has achieved the designed milestones, including (a) 
developing an online automatic RWIS validation tool, (b) validating and refining the icy road 
forecast algorithms, and (c) developing an in-depth understanding of the formation of road ice. 
Specifically: 


The overall hourly road ice forecast accuracy improved from an averaged 62% to 
82%. This statistical result is based on the comparisons between the hourly 
forecast for all days in a month with the corresponding RWIS observations across 
Montana. 


Though progress on road ice formation mechanisms and forecasting capabilities has been 
made, research on black ice was very limited due to a lack of black ice identification and funding 


' Product supported through research funding. SpringGem may not have resources to actively maintain webpages 
listed here and in the rest of this report upon project completion. 


availability. Future field experiments could help collect information on conditions leading to 
black ice occurrence on US roadways. More development on model forecast schemes for black 
ice formation is also needed. 


RWIS site observations collection 


Sub-hourly observations for January 1 — March 31, 2020 and October 1, 2020 — March 
31, 2021 were collected from 72 RWIS sites through the RWIS data access webpage 
(www.rwis.mdt.mt.gov). Data was cleaned by removing errors and/or missing values. Data from 
one RWIS site was omitted due to sensor failure. Removed error data and flagged data when the 
sensor status reported “chemically wet”. This study also developed a JavaScript code to acquire 
data from the RWIS sites and process that data automatically. 


RWS observation storage in database 


RWIS data was saved in an online cloud-based database, which made it possible to 
conduct web-based automatic data analysis for all RWIS sites. The online database also helped 
others, in particular MDT managers, to check the performance of the IcyRoad forecast by 
comparing with RWIS sites for specific regions of interest. This step requires coding in Python, 
JavaScript, and database management. The webpage is available at: 
https://sgweather.com/INTERNAL/IcyRoadValidation/. 


Data Analysis 


New understandings of road ice formation conditions resulting from this project are 
shown in Figure 1| and include the following four points: 
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Figure 1: RWIS instrumentation observations on Bozeman, MT, January 2020. The x-axis is the data entry from Bozeman RWIS 
site, for ice-occurring times only. RWIS data has interval of about 5 minutes. Data is available at https://rwis.mdt.mt.gov. 
Bozeman Pass has three sensors and this study used the WB lane MP 321.3 since the other two have errors. 


1) Urban regions had black ice* much more frequently than other regions (results from 
other sites not shown). Bozeman, MT for example, is a city in southwest Montana with a 
population of approximately 49,631 in 2019. The urban heat island effect led to no 
rainfall/snowfall in the city during January and March 2020, but ice on road surfaces still 
occurred (Figure 1) ~80% (24 days versus 31 days) of the time, due to strong radiative cooling at 
the surface, low surface temperature, and abundant water vapor in the atmosphere. This finding 
suggests that a different city-road ice scheme from other land-cover-surrounding road regions, 
specifically rural or mountain regions, is needed. 

2) Road surfaces were generally warmer than the 2-meter air temperature by 1-10 °F (- 
17.2 - -12.2 °C) but had clear diurnal variations, partly due to anthropogenic heat flux from 
traffic and the modified physical properties of the paved materials (e.g., albedo, thermal 
emittance, thermal conductivity, specific heat, and surface roughness). This is in contrast to the 
previously well-known understanding that at night the land surface is cooler than the 2-meter air 
temperature due to radiative cooling at the ground. It suggests, again, that roadways are a unique 
land surface category and should be simulated in weather forecast models using a different 
physical scheme to represent their radiative and thermodynamic features. 

(3) Soil temperature varied in a much smaller magnitude than road surface temperature” 
and 2-meter air temperature, in both diurnal range and monthly scale. This may be due to the 
large heat capacity and hydro-conductivity of soil layers underlying RWIS sites. Such a small 
change in soil temperature indicate that this variable cannot reflect the rapid changes of the 
ground surface energy modification or be used to forecast road ice formation. 

(4) Road surface temperature and 2-meter air temperature had a high correlation coefficient up to 
0.83 in January 2020 for Bozeman Pass. Air temperature varied more rapidly than road surface 
temperature, for example during the 1600-2000 period of time (x-axis, Figure 2), due partly to 
advection. These results imply that deriving road surface temperature from air temperature 
together with other weather information is possible but may miss detailed changes. A previous 
study (Jin et al. 1997) discussed that during cloudy days, road surface and air temperature are 
close to each other, but that during sunny days, road surface temperature is higher than air 
temperature. For highways, a human-induced environment, road surface temperatures were 
higher than air temperature for most day and night times but could be lower when synoptic 
advection warm air occurred. As a result, the road surface energy budget differed from any other 
land surfaces and balanced among downward longwave radiation, upward surface-emitted 
longwave radiation, upward sensible heat flux from surface to atmosphere (opposite direction 
from normal vegetative/bare-soil land surfaces), and upward ground heat flux (e.g., from soil 
layers to surface), assuming the latent heat flux is negligible for the water-proof roads. 
Therefore, this study developed a new road energy scheme which better simulated road 
temperature in current WRF/land-surface model forecasts. One critical finding is that when road 


There are two kinds of land surface temperatures: 2-meter air temperature and land surface skin temperature. Two- 
meter air temperature is the one traditionally used in the transportation sector and is measured by a weather station 
with a thermometer in a Stephen shelter built 1.5 - 2 meters above the ground. Skin temperature, a radiometric 
temperature by definition, is measured via remote sensing from thermal infrared emission as a function of ground 
temperature. The road surface temperature measured by RWIS is close to skin temperature. Road surface 
temperature is not simulated and thus not forecasted in the weather model (Jin et al. 1997, Jin & Dickinson 1999, 
2000; Jin et al. 2005; Jin & Shepherd 2008; Jin 2012; Jin et al. 2014). Skin temperature is the one variable that 
actually determines ice formation on the road. Skin and 2-meter air temperatures have a high correlation due to 
surface boundary-layer processes but significant differences in terms of magnitude and physical meaning (Jin 2000). 
Road surface tempera-ture observed in RWIS is equivalent to road skin temperature. But to be consistent with other 
documents, this report uses road surface temperature instead of road skin temperature. 


surfaces have ice, road surface temperatures were equal to or below 32 °F (0 °C)? (Figure 2). Air 
temperatures could concurrently be higher than 32 °F (0 °C). This finding may suggest that if 
road surface temperatures could be forecasted accurately, it would improve an ice road forecast. 
This finding may also correct a common misunderstanding among the public: often when reports 
of a road ice-induced vehicle accidents are made, they suggest that road icing occurred even 
when temperatures were above freezing. In fact, although 2-meter air temperatures may be 
beyond freezing, the road surface temperatures were always equal to or below 32 °F (0 °C), as 
shown in Figure 2. 


Bozeman Pass, Bozeman, MT, Icy Road, January 202 
60 7 


40 


20 


Road surface Temperature (F) 


-20 10) 20 40 60 
2-m Air Temperature (F) 


Figure 2: The road surface temperatures and 2-meter air temperatures for Bozeman Pass for January 2020. Only data indicated 
as “icy status” by RWIS sensor were analyzed in this plot. The pink line is where skin temperature equals to air temperature. The 
dashed lines parallel to x-axis and y-axis were 32 °F (0 °C), respectively. 


Road surface temperatures (Troaa) were higher than the 2-m air temperatures (Tair) most of 
the time, since most points in Figure 2 occur above the pink line, where Troad and Tair are equal. 
However, air temperatures higher than Troad also occurred. This is a major challenge in 
forecasting: a better understanding is still needed when air temperatures can be higher or lower 
than road surface temperatures and hopefully form an approach to predict such a Troad/Tair 
scheme. This scheme would be related to local surface-layer turbulence, convection, and 
moisture evaporation. 


3It is unclear to the researchers how exactly RWIS determines road status (ice warning, wet, chemically wet, etc). 
RWIS defines “wet” as “Continuous film of moisture on the pavement sensor with a surface temperature above 
freezing (32 °F or 0°C).: see 

https://rwis.mdt.mt. gov/scanweb/SW Frame.asp?Pageid=SiteOverlayMap& Units=English&Groupid=267000&Siteid 
=267007 &Senid=&DisplayClass=Java&SenType=All&Mapid=571. This pre-set, instrument-level threshold may be 
a reason for no ice found above 32 °F. Research is needed to measure road surface temperature before, during, and 
after the ice occurrence and when the road is dry. There may be uncertainty in RWIS data and thus results based on 
this data. 


The MacDonald Pass site in Jefferson County is a high mountain region with an elevation 
above 6,320 feet (1,926 m) The road at this site was icy for >90% of days in January 2020 
(Figure 3). 
Snowfall was also monitored in most of the days on this site (77%), partly due to convection 
induced by orographic lifting and partly due to wind-blown snow from nearby regions being 
measured by the RWIS sensor. Furthermore, the north-facing aspect of the mountain received 
relatively less irradiance resulting in snow persisting as a water source for road ice. Road 


MacDonald Pass, Butte, MT Icy Road Temperatures 


100 = 100 
¥. 
80 80 
60 60 


40 


Py 
=) 
Relative Humidity (%) 


Temperature (F) 


20 


OFT sata 2-meter Air Temperature ; i corr. Air and Road Temperatures = 0.85 9 
x Precipitation occur with precipitation case = 4413, (77.5%) 
—------ dew point 
RH 
-20 -20 
O 2000 4000 6000 


Obs Records in January, 2020 


Figure 3: The same as Figure 1 except for MacDonald Pass site, Butte. The road surface temperature, 2-meter air temperature, 
dew point, relative humidity, and precipitation presence are presented. Dashed line is for 32 °F (0 °C). 


surface temperatures were highly correlated with 2-meter air temperatures with a correlation 
coefficient of 0.85. For this mountain site though, road surface temperatures were often lower 
than 2-meter temperatures, due partly to the 90% of time when ice was present on the ground. In 
addition, road surface temperatures were all below 32 °F (0 °C) for ice-occurring times, which is 
represented by the dashed line in Figure 3. 

Comparing with the Bozeman Pass site (Figure 1), the different Troad/Tair relation suggests 
that local calibration is needed in order to forecast Troaa correctly from Tair for a given region, since 
local elevation, cloud cover, and geographic conditions directly affect Troad/Tair interactions. Land 
surface has high heterogeneities. The degree to which one observation of roadways could represent 
the local conditions is a question still needing to be studied. 

Relative humidity (RH) varied from 23-100% at MacDonald Pass when ice occurred. 
This indicated that in this region, ice formed on the road regardless of whether the RH was high 
or low. Therefore RH, or dew point which could be derived from RH and 2-meter air 
temperature, was not a strong indicator for road ice occurrence. 


Data Analysis Key Conclusions 


(a) Synoptic weather processes determined the overall precipitation, air temperature, and 
wind conditions at each site. Local land cover, geographic features, and land-water 
fractions also contributed to road surface and air temperatures. 

(b) Black ice can form within a range of relative humidity from 20-100%. 

(c) The road surface temperature was the key parameter for road ice formation. 

(d) No ice was found when the road surface was above 32°F (0 °C). This 32 °F (0 °C) 
threshold for road surface temperature should be a strong indicator for no occurrence of 
road ice. Therefore, predicting road surface temperature with adequate accuracy is 
critical. 


SECTION I 


Online Validation Tool* Development for Data Process Efficiency 


Internet-of-Things technologies makes it possible to rapidly access, visualize, and 
analyze road ice forecasts. Unlike traditional methods for downloading and analyzing data site 
by site from a local computer, the online, automatic validation tool can access 72 RWIS sites 
across Montana and analyze the forecast and measurements in a few seconds (Figure 4). 

In addition to speeding up validation analyses, the tool's automatic, hourly validation can 
also be used to monitor the performance of both RWIS observations and forecast models. See the 
daily updated tool at http://sg¢-weather.com/INTERNAL/validation/ (Figure 4 is an example, 
comparing RWIS hourly measurement with model forecast). 

The upper panel of Figure 4 allows users to select which RWIS site to examine. Once a 
site is selected from the dropdown menu, the location of the site is highlighted in red on the map. 
The hourly surface temperature and road ice forecasts are concurrently accessed, compared, and 
visualized in the lower panels, by comparing with the mean RWIS hourly observations. 
Furthermore, subsets of the data can also be analyzed using a built-in crop function in the tool. 
This tool serves to enhance transparency of validation. 

Hourly, daily, and monthly statistical analyses were conducted by comparing the 
National Oceanic and Atmospheric Administration (NOAA) North American Mesoscale 
Forecast System (NAM) model-based IcyRoad forecast with the MDT 72 RWIS sites. 
Specifically, two variables were analyzed: ground surface 2-meter air temperature and ice road 
status. The forecasted 2-meter surface air temperature had high accuracy with a correlation 
coefficient higher than 0.80 for >70% of sites. The ice road status was relatively lower, at 0.64 
for all sites. The ice road status was lower partly because the algorithm (IcyRoad2) needs further 
improvement on the black ice formation scheme and partly because the MDT highway 
maintenance program applied anti-icing material (i.e., “chemically wet”) in an effort to help 
prevent icy roads. Therefore, this forecast differed from the RWIS observations on ice status. 


4 Separate online tools were designed for each IcyRoad forecast algorithm, which helped to assess accuracy of each 
algorithm. It is recommended that MDT use https://sg-weather.com/INTERNAL/IcyRoadValidation/, which not 
only compares RWIS hourly data with model (IcyRoad3 algorithm) forecast, but also has statistical data analyses for 
day versus night time, ice and no ice. 
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Figure 4: Online automatic validation of IcyRoad2 forecast with MDT RWIS sites. (a) shows the 72 RWIS sites with the selected 
site in red. (b) Middle panel is the IcyRoad2 forecasted 2-m air temperature (in green) comparing with the RWIS sensor observed 
road surface temperature (in red). (c) The lowest panel is the IcyRoad2 forecasted icy status and RWIS report-ed road status. 
The website is available at http://sg-weather.com/INTERNAL/validation/. IcyRoad2 is one of forecast algorithms developed and 
examined. 


Note that in Figure 4, the NOAA NAM218 weather forecast was used as the base 
forecast. NAM218 has a 12-km spatial resolution and hourly temporal resolution and is available 
for free to the public. In general, NAM218-IcyRoad2 forecasted the diurnal variation of surface 
air temperature well, but generally underestimated the actual temperature by about 1-5 °F (-17.2 
to -15 °C) during January—March, 2020. Namely, the red line was almost consistently above the 
green line. Consequently, the ice road forecast, which was determined by the surface air 
temperature and road surface temperature, was affected with an accuracy of about 0.67 for this 
specific site (Ninemile I-90 MP). Nevertheless, different sites had different road ice forecast 
accuracy. 


This online tool can also efficiently disclose the relationship between the RWIS observed 
road surface temperature and road temperature (Figure 5). First, both road surface and air 
temperatures have clear diurnal variation, following the diurnal variation of surface insolation. In 
addition, at this station and for this period of time, road surface temperatures were almost always 
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Figure 5: Online tool for road surface temperature and air temperature relationship. Data shows hourly RWIS observations for 
the Ninemile I-90 MP site from February 25-March 13, 2020. 


higher than the 2- meter air temperature. For clear days, surface temperatures could be higher 
than air temperature by about 5 °F (-15 °C) (March 3, 2020) to 30 °F (-1.1 °C) (February 25, 
2020). Nevertheless, at night, road surface temperatures were much closer to the air temperature, 
but still slightly higher. Even when the road was icy, surface temperatures were still higher than 
the air temperature. For example, on February 27, 2020 from 0200 to 1900, the RWIS site 
indicated ice with the road surface and air temperatures were very close to each other. 


Evaluation of Model Forecast for Surface Air Temperature 


Local geographic conditions evidently affect 2-meter air temperature (Figure 6). 
NAM218 was able to simulate these effects well. In the Missoula region, 5 RWIS sites were in 
service: Lookout Pass, Bearmouth, Trout Creek, Greenough Hill and Ninemile. Each RWIS site 
had one or more sensors. It was examined how weather forecasts performed for these different 
sites. First, due to the relatively close proximities of the sites, they had very similar diurnal 
variations. This suggested similar sky conditions, clear or cloudy, rainfall or no rainfall for these 
sites. Second, different sites had different diurnal ranges of air temperature. For example, around 
step 50 (x-axis), Bearmouth had a much larger diurnal range - 268 K (22.7 °F, -5.2 °C) at noon to 
254 K (-2.5 °F, -19.2 °C) at the lowest time around sunrise next day. In addition, the Ninemile 
site almost always had higher temperatures than other sites. Further analysis showed that soil 
moisture, land surface albedo, and vegetation fraction near the site all contributed to a different 
diurnal temperature range. These results suggested that the model can capture the local 
geophysical effect on air temperature reasonably well. 

For the Great Falls region (Figure 7), similar diurnal variations occurred at each RWIS 
sites. More importantly, Bowmans, a site with an elevation of 4,301 feet (1311 m), always had 
the lowest temperatures at night. For these sites, the sky conditions were not uniform. For 
example, at 290" hour (x-axis), Roger’s Pass had the coldest two days without the temperature 
increasing during the day, while other sites had a small diurnal range, indicating more severe 
cloud conditions with precipitation. 


10 


NAM 218 forecast, February, 2020 
290 T T T T ———— : T T T 
— Lookout Pass, Missoula 
—— Bearmouth, Missoula 
~~ Trout Creek 


Greenough Hill, Missoula 


280 , _ Ninemile, Missoula P 
8 M 7 | nee \ 
8 y | J | j \ \ 
zz \ Nw A } A f\ Wl h / A a \ \ 
| ASW ep AT Ww hp 4g | ¥ e } 
= 270 fh Aad Wi) ea hk | Ve Mi) 
= \ ( ” y s 
: MTR RAINY NS 
2 S 4 Vi bh \ \ \ \ 
& ' : W | \ 
5 Nu N \ 

260 y \ 

2x0 : l ! l — ' a ! 

0 100 200 300 400 500 600 700 


Hours of days in February 
Figure 6: Weather forecast of 2-meter air temperature variations for the Missoula RWIS sites locations. X-axis is hours from 
February 1, 2020-February 28, 2020. Y-axis is air temperature in unit K. 


Surface temperature heterogeneity was the result of both local geographic conditions as well as 
cloud fraction and rainfall heterogeneities. 
NAM 218 forecast, February, 2020 
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Figure 7: Same as Figure 6 except for RWIS sites at Great Falls. 


Evaluating IcyRoad Forecasting Algorithm 


All weather models have uncertainties given that models, by design, reflect the current 
understanding within the scientific community on physical and dynamic weather processes. The 
current understanding of road ice formation is far from complete. To gain insights on how to 
improve road ice forecast accuracy, a few numerical algorithms were designed, depending on 
various turbulence theories and data assimilation methods. 

Results below (Figure 8) were based on the weather forecast produced by NOAA's NAM 
model for January-March 2020. NAM is one of the weather forecast products used in this 
project, with other weather models including commercial Openweather.org, NOAA's Global 
Forecast System model with spatial resolution of 15-km, and SpringGem’s developed weather 
forecast. Surface temperature was found to be one of the most important factors for road ice to 
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form. The temperature threshold is therefore critical to road ice forecasting. Designed 
temperature thresholds of 32 °F (O °C), 28 °F (-2.2 °C), and 24 °F (-4.4 °C) in the algorithm and 
ran the forecast model to calculate its accuracy for different temperature thresholds. 


Nam 218 Ice Status Forecast Algorithms, January 2020 
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Figure 8: NAM-IcyRoad2 accuracy for road ice status forecast in January 2020 at RWIS sites. X-axis is the RWIS site names and 
y-axis is total accuracy percentage for ice status forecast. 32, 28, and 24 represent surface temperature thresholds at 32°F (0° 
C), 28°F (-2.2 °C), and 24°F (-4.4 °C), respectively. 


Total ice forecast accuracy percentage, the y-axis, was calculated as follows: for each 
hour from January 1-January 31, 2020, the NAM-IcyRoad2 forecasted road ice status was 
compared with RWIS hourly observation for that specific site. Agreement between the two was 
counted as "1", and disagreement counted as "0". The sum of "1" divided by the total forecast 
hours is the total ice accuracy percentage. Therefore, in Figure 8, the highest accuracy was 94% 
for Homestake Pass I-90 associated with a temperature threshold of 32 °F (0 °C). 

A valuable understanding of road ice forecast accuracy was gained (Figure 8). Different 
sites had different sensitivity to temperature thresholds. This so-called local calibration are 
defined by the geophysical conditions of the site. For example, the sites in the Missoula region 
(site numbers 150000-150005), such as Trout Creek had the highest accuracy at the 24 °F (-4.4 
°C) threshold (>80%) and lowest for the 32 °F (0 °C) threshold (<60%). The nearby Greenough 
Hill site however, was most accurate at the 32 °F (0 °C) and 28 °F (-1.1 °C) thresholds (~76%) 
but lowest in the 24 °F (-4.4 °C) case (< 60%). Furthermore, the 32 °F (0 °C) threshold was the 
most acceptable for almost half of the sites (27 sites vs. 57 sites analyzed). Sites in the Miles City 
(563000, Government Hill, Ingomar US12), Bozeman (564000, Bozeman Pass I-90; Norris Hill 
US287), and Havre regions (629000, Inverness US2 MP 337.6; Loma US87 MP 53.1; Sweet 
Grass I-15) were most accurate for the 24 °F (-4.4 °C) threshold. 

In order to determine if the information detailed above varies monthly or seasonally, 
another month of data (March 2020) was further analyzed (Figure 9). 
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Nam 218 Ice Status Forecast Algorithms, March 2020 
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Figure 9: Same as Figure 8 except for March 2020. 


Obvious consistencies were observed for March 2020 road ice forecast accuracy (Figure 
9). First, 32 °F (0 °C) and 24 °F (-4.4 °C) were the two most accurate thresholds. The Missoula, 
Miles City, Bozeman, and Havre sites were still most accurate at 24 °F (-4.4 °C) and other sites 
from Butte (267000, i.e., MacDonald Pass) to Lewistown (268000, 1.e., East of Denton MT81) 
were most accurate at 32 °F (0 °C) (except for the Garrison site in Butte). The temperature-only 
algorithm had an average accuracy of 62%, with some sites close to 80% but some sites below 
60%. This suggests that a temperature-only forecast algorithm is not adequate. This study need to 
use other information, such as previous precipitation data, clouds, and local geographic 
conditions (1.e., elevation, vegetation fraction, distance to water source). 

We also need to improve surface temperature forecasts. The road ice formed under a 
combination of road surface temperature, road air temperature, wind, and water availability. 
Obtaining the correct road surface temperature from 2-meter air temperature remains a 
challenge. Another challenge is in identifying from the road surface temperature the threshold at 
which ice can form on the road surface. 


Further Research on Ice and Black Ice Formation 


The Bearmouth site (ID 150002) is a mountain-like site with an elevation of 3,901 feet 
(1,189 m) and close proximity to the road. For this site in January 2020, road ice occurrences 
with current precipitation cases totaled 803, or 16% of total ice occurrence cases. Nevertheless, 
when previous precipitation is considered, 55% of road ice occurred when precipitation occurred 
in the previous 24 hours (Figure 10). 
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Figure 10: RWIS site Bearmouth, Missoula for January 2020. 


The road surface temperature and 2-meter air temperature correlation coefficient was 
0.92. Air temperature can be much higher than road surface temperature during the daytime, 
partly due to the fact that ice-covered road surfaces do not heat up quickly. Nevertheless, air 
temperature could also be significantly lower than the road surface temperature (>10 °F (-12.2 
°C) difference). Furthermore, it seems that on clear nights, air temperatures were lower than the 
road and black ice could form relatively easily. 

When roads were icy, the surface temperatures were always below 32 °F (0 °C) (Figure 
11). Nevertheless, 2-meter surface air temperature could be above 32°F (0 °C). Therefore, 32°F 
(0 °C) surface temperature should be a critical threshold for road ice not occurring. Again, when 
there was precipitation for icy road conditions, the surface and air temperatures were close to a 
linear relation (1.e., close to the pink line which indicates the two temperatures were equal). 
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Bearmouth Missoula, MT, Icy Road, January 2020 
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Figure 11: Hourly RWIS instrumentation observations for Bearmouth site of 2-meter air temperature and road skin temperature 
for icy road occurrence times, from January 1—January 31, 2020. The green markers were the cases of precipitation at the time. 


Elevation Effects on IcyRoad Forecast Accuracy 


Elevation affects the ice forecast accuracy of ice formation on roadways. Figure 12 shows 
a negative correlation coefficient of -0.25 between RWIS site elevation and icy status forecast for 
the site. Although this negative coefficient may not be statistically significant, elevation is still one 


of the secondary factors. 
January, 2020 IcyRoad Forecast on Montana RWIS Sites 


—Icy Forecast Accuracy 
— Site Elevation (unit:feet) 


pee 


Figure 12: Elevation (unit: feet, blue line) for each RWIS site and Ice Accuracy Percentage (red line) for 
January, 2020 at temperature threshold t=24°F (-4.4 °C). IcyRoad2 algorithm was examined in this Figure. 
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Surface Temperature Accuracy Effects on IcyRoad Accuracy 


Figure 13 shows that overall the NAM-forecasted 2-meter air temperature has a 
correlation coefficient for most of the sites above 0.70, with 45% above 0.80. Nevertheless, the 
icy road accuracy (y-axis) was not linearly related to temperature forecast (x-axis). But the more 
accurate the temperature forecast, the better chance to forecast road ice correctly. This suggests 
that an improved accuracy of surface temperature forecasting is needed. Road ice is a function of 
road surface temperature, road air temperature, wind, and water availability. 


January 2020, Montana RWIS Validation on NAM forecast 
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Figure 13: Correlation of forecasted air temperature and ice forecast accuracy. NAM-IcyRoad2 was evaluated in this plot. 


SECTION II 
RWIS validation on Refined IcyRoad3 Forecast 


Assimilating the past 6-12 hours of ground observations into weather forecasts efficiently 
improved the IcyRoad forecast. IcyRoad3 was designed to assimilate RWIS instrumentation 
observations of air temperature and rainfall information. Hourly measurements from RWIS sites 
across Montana were used to compare with the corresponding IcyRoad3 forecast (Figure 14). 
IcyRoad Accuracy Percentage, the y-axis in Figure 14, was calculated by (1) if for a specific 
hour, the IcyRoad3 forecast was the same as the RWIS ice status measurement (icy, not icy) then 
it was counted as "1", otherwise, it was counted as "0"; (2) during all 24-hours, 28 days in 
February, summed all "1" cases; and (3) divided by all forecast entries (1.e., 24 hours multiplied 
by 28 days), to get the Icy Status Forecast Accuracy Percentage. Based on this definition, in 
February 2021 (Figure 14): 


e The accuracy for ice status was 62-96%, with the all-site-averaged accuracy being 
82%. An averaged 82% accuracy is very promising. 

e There are only three sites with an accuracy below 70%. The lowest two sites were 
Bozeman (62.56%, site ID 564000) and Gary Cooper Bridge (67.83%, site ID 
301004.1). One is an urban site and the other a bridge. Both sites had errors 
(Figure 14) with RWIS being icy but IcyRoad forecasting “no ice”. The 3" site is 
Monida Pass (69.76%, site ID 267001.1). This is a mountain site with an 
elevation of 6,824 feet (2,080 m). 
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e Almost all sites underestimated road ice, namely, that RWIS observed ice but 
IcyRoad did not forecast it. This means there was systemic uncertainty in the 
forecast model that should be removed. 


Therefore, to improve the accuracy of the forecast, this study need to: 


e Perform in-depth analysis on ice formation for the urban, bridge, and mountain 
sites. This study may need to modify the forecast algorithms for these categories, 
and 


e Refine the thresholds of this ice detection algorithm. 


Icy Accuracy Percentage, February 2021, IcyRoad3 Algorithm 
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Figure 14: Icy Status Forecast Accuracy Percentage across all Road Weather Information Sites (RWIS) for February, 2021, X- 
axis is the site names. Y-axis is the ice status forecast accuracy percentage; calculation method is described in the text. 


Further Ice Formation Study on Bridges using Online Data Analysis Tool 


In order to rapidly access, analyze, and visualize RWIS sites for validation of the 
IcyRoad3 forecast, an online validation and statistics analysis tool was developed at https://sg- 
weather.com/INTERNAL/IcyRoadValidation/. This tool combines Internet-of-Things, an online 
database, cloud computing, and statistics analysis into an advanced research tool. 

As an example, this tool was used to understand the Yellowstone Bridge (site ID 
263000.0) road surface conditions for ice versus no ice and day versus night, respectively. In icy 
cases, road surface temperatures at the Yellowstone Bridge site were always equal to or below 32 
°F (0 °C) (Figure 15). This means that no ice can be formed when the road surface temperature is 
above 32 °F (0° C). Since all weather stations only measure 2-meter air temperature and weather 
forecasting models only provide 2-meter air temperature information, it is critical to know the 
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road surface temperature. The least-square linear relationship between Troad and Tair is calculated 
with the following equation: 


Tréad = 0.95 Tair + 7.21 Eq. (1) 


Yellowstone River Bridge |-90 MP 452.3 - App 


; 


Icy 
= Wot Icy 
— y = 0.956 + 7.24 
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Road surface Temperature (F) 


—20 oO 20 40 6a oO 
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Figure 15: The 2-meter air temperature and road skin temperature relationships for November 1, 2020 — to March 31, 2021, for 
Yellowstone Bridge (263000.0). The blue dot represent no-ice time, and the red dot represents for icy. 


Troad-Tair relation, Eq (1), is only a statistics-based relationship and different sites have 
different equations due to local land cover and micro-meteorological conditions. Here only one 
year (January-March, 2021) of data was used in the analysis. If a 10-year dataset were to be 
analyzed, the statistical analyses would be statistically significant. 

The Troad-Tair relation has been a challenge that research communities tried to 
address in the last two decades but there is still a lack of deep understanding. This relation is 
determined by surface layer turbulence which is very difficult both to measure and to simulate. 
Using the online tool, can the Troad-Tair relation can be obtained, as a first-order estimation. 

The difference between skin and air temperatures is determined by solar radiation 
(longwave radiation emitted from the atmosphere to the ground and from the ground back to the 
atmosphere), sensible and latent heat fluxes, and ground heat flux including anthropogenic flux 
due to traffic. Therefore, this relation has a diurnal cycle and thus may need to be split into 
daytime and nighttime components, respectively. 

During the daytime, the relation is Troad =0.97Tairt+7.52, and at night, Tskin = 0.95Tairt+6.97 
(Figure 16). This suggests that at night, with a 1 °F temperature decrease on Tair, Troad reduces by 
0.95°F, slower than Tair. Furthermore, when Tair was lower than 14 °F (-10 °C), the road always 
had ice. This value could be used later as an important threshold in the road ice forecast 
algorithm. 
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Figure 16: Same as Figure 15, except for (a) daytime (6 AM — 6 PM) and (b) nighttime (7 PM to 5 AM), respectively. 


Ice, including black ice, forms most frequently at night or during the early morning due 
partly to road surface radiative cooling. Therefore, further study into what conditions lead to 
nighttime ice formation is needed. Again, all road ice formed only when the road skin 
temperature was equal to or below 32°F (0 °C) (Figure 17a). 

Using the relation, Troad =1.02Tair+ 12.86, during the night at Yellowstone Bridge, Tskin 
dropped faster than Tair. In contrast, when roads had no ice (Troad =0.98Tair+0.91), Tskin dropped 
slower than Tair (Figure 17b). The faster Troaa decrease may be a reason for road ice formation 
but may also be a consequence of road ice. Nevertheless, once Tair dropped below 14°F -10 °C), 
road ice always formed. 

From November 1, 2020 to March 31, 2021, total data entries for the Yellowstone Bridge 
site was 2,108, with a total of "ice" reports being 473 and total "no ice" reports of 1,635. Among 
this data, the correlation coefficient between 2-m air temperature relative to road temperature for 
no ice cases was 0.78; the correlation coefficient for total 2-m air temperature relative to road 
temperature was 0.93 for ice cases. 
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Figure 17: Same as Figure 15 except for (a - top) ice cases at night and (b - bottom) no ice cases at night, respectively. 
Temperature unit °F. 


Similar research has been conducted on other bridges but are not presented here as those results 
are easily obtained via the online analysis tool. 


Refined IcyRoad3 Forecast in March 2021 — Current Status 


After refining the forecast algorithm by updating temperature thresholds and Troad-Tair 
relations for urban, bridge, and mountain sites, the IcyRoad3 algorithm evidently improved 
forecast accuracy (Figure 18). All sites except one, Swan Lake South (accuracy 61.03%, site ID 
628005.0), had an accuracy above 80%, with the statewide average accuracy as high as 88%. 
Furthermore, the systemic error disclosed in Figure 17 was corrected, with part of the error due 
to RWIS "Ice" forecasted as "No Ice" and part due to RWIS "No Ice" but forecasted as "Ice" 
(green bar and red bar are now close to each other) 
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Figure 18: Same as Figure 14, except for March, 2021. 


Conclusions 


The Statewide Icy Road Forecasts for Montana was developed and is presented at 
https://sg-weather.com/INTERNAL/Icyroad3/ (Figure 19). The road ice forecast, RWIS site 
observations, and weather stations are presented on the same page to enable cross-checking of 
information. The forecast lead time is currently 24-hours but can be extended to 72 hours if 
desired by the needs of MDT. 


General conclusions for Task | are as follows: 


a. Monthly hour IcyRoad3 forecast accuracy is above 80% (see Figure 17 and 18). For 
certain regions and/or under certain weather conditions, the accuracy may be higher or 
lower. 

IcyRoad3 covers most highways in Montana with reasonable accuracy. 

c. The longer the forecast lead time, the lower accuracy of IcyRoad forecast. This is a 
common limitation for all weather forecasting, since weather forecasts are sensitive to the 
initial conditions and thus involve inherent limits on predictability. Given this, re-checks 
(i.e., 3-hourly) of the IcyRoad3 forecasts are recommended. 

d. The results suggest that a different city-road ice scheme from other land-cover- 
surrounding road regions, specifically rural or mountain regions, is needed. 

e. Obtaining the correct road surface skin temperature from 2-meter air temperature remains 
a challenge. Another challenge is in identifying the threshold at which ice can form on 
the road surface given road surface skin temperature. 
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IcyRoad Forecast for Montana 
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Figure 19: Sample of IcyRoad3 for state Montana highway ice status forecast for February 5, 2021. The blue marker is icy, the 


red marker is no ice, the black marker is RWIS sites location. Clicking RWIS sites shows observations. The small, pink pins are 
weather station locations and observations. 
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Task 2 Validation of the Ice Formula 
SECTION I 


Introduction 


Black ice mechanisms are one of the most meteorologically challenging processes to 
model and more observations are needed to understand black ice forming meteorological 
conditions, locations, and timing of events in order to forecast it well. Since the road ice 
observations, RWIS sites, are located only in limited areas along roads, more observations 
between RWIS sites are needed. Studies have been done using additional sensors such as near 
infrared cameras (0.9 - 2.3m) with evidence of possible black ice detection although wet and 
black ice surfaces had similar sensor responses and new technologies were recommended to 
improve detection accuracy (Jonsson et al. 2015). New technologies are definitely needed in this 
project to validate IcyRoad forecasts since visible cameras and RWIS data cannot discern the 
difference between wet and black ice surfaces. Therefore, the team has been developing a drone- 
based remote sensing technology to detect black ice via hyperspectral camera launched on an 
unattended aerial vehicle (UAV). This detection approach has been intended to help us to 
improve physical modeling of black ice which has been ignored in the existing Weather 
Research and Forecasting model (WRF/land-surface) and could validate the black ice forecast 
component in the Icy Road Forecast and Alert system. There is evidence of hyperspectral 
cameras with visible and near infrared capabilities (0.4 — 1 um) being used to identify road 
defects (Abdellatif et al. 2019) but no one has yet attempted to distinguish road surface 
conditions using hyperspectral cameras. This project has used a hyperspectral camera capable of 
distinguishing 281 spectral bands at 2.1 nm resolution. As a result, an increase of accuracy for 
IcyRoad forecasts has been expected. 

Black ice is hard to detect from the visible bands (i.e., 0.4-0.7 um, the wavelength 
of solar radiation spectrum that human eyes use for sight, see Figure 20). Nevertheless, 
an ice sheet, thin or thick, has spectral reflectance (ra) significantly different at near 
infrared bands. By using the combination of visible bands and near-infrared bands (NIR), 

a thin ice-covered surface may be differentiated from asphalt surfaces (1.e., road 
surfaces). Specifically, in which RED is the radiance for the red band. 


Road Ice Index 6 = (NIR-RED)/(NIR+RED) (1) 
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Figure 20: Spectral reflectance curve for asphalt, water, snow and vegetation. 


The Red band is the edge of the visible band at close to 0.620 — 0.670 um (e.g., MODIS 
band 1, King et al. 2004). NIR is the near-infrared region of the electromagnetic spectrum 
(from 0.78 tum to 2.5 um) and this analysis started using 0.841 — 0.876 um (e.g., MODIS 
band 2) to follow NASA MODIS instrument spectral bands. Theoretically, for snow and 
ice, ra reduces from RED to NIR and thus 0 is negative; while for asphalt surfaces, ra 
increases from RED to NIR and thus has a positive 8. As a result, icy roads can 
theoretically be detected using a combination of Red and NIR bands. Nevertheless, this 
all depends on the accuracy of the hyperspectral calibration and the need to use an 
additional camera capable of longer bands in the NIR was explored, such as MODIS 
band 6 (1.628-1.652 um) since research showed evident spectral difference for dry, wet, 
snow, and icy roads at this band (Hall et al. 2002; Jonsson et. al 2015). 


SECTION II 


Data Collection 


Main objectives for the data collection were; 

a) First, the team collected baseline observations of spectral 
signatures on dry, wet, ice, and possibly black ice on asphalt in a 
controlled environment. The Montana State University Sub-zero 
laboratory was used to control temperature and wind conditions 
while setting up various scenarios on asphalt samples (i.e., snow- 
ice, old-snow-re-freeze, and black ice cases). 

b) Second, the laboratory measurements were repeated (using a tripod 
mount for the camera) at the Missoula International Airport on 
runway 26/8 near taxiway C during times when the runway was 
shut down in the winter. This gave the team a secure site with 
access to asphalt for conducting imaging under a variety of 
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conditions (dry, wet, and ice) while avoiding roadway/runway 
traffic. 

Cc) The same hyperspectral camera from the subzero laboratory was 
used on an uncrewed aerial vehicle to detect road ice conditions at 
varying altitudes when no traffic is present on UM property. 


Experimental Setup 


September 21‘ - 25", 2020 over 70 GB of data was collected in the Montana State 
University sub-zero laboratory using two different hyperspectral cameras, UM's Pika L and a 
borrowed Resonon's Pika L 320 whose range is farther into the NIR than the UM Pika L. This 
data resides in UM’s box cloud storage. Each spatial point (pixel) in the hyperspectral images 
represents a continuous curve of incoming light intensity versus wavelength. The data can also 
be interpreted as a stack of images, with each layer in the stack representing the scene at a 
different wavelength, this “stack” of two-dimensional images is referred to as a “datacube” of 
three dimensions. 

Each day a camera setup on a tripod mount was used to take imagery at air temperatures 
from 24 °F (-4.4 °C) to 40 °F (4.4 °C) at 2-degree increments. These measurements were 
repeated from 40 °F (4.4 °C) to 24 °F (-4.4 °C) at 2-degree increments. See Figure 21 for 
experimental setup. Once an image was recorded Resonon’s software was used to calibrate the 
data based on the calibration panel. Both the “raw” and “calibrated” data cubes were saved. Then 
ArcGIS was used to pull out individual wavelength matrices for ingestion into the ice index 
formula. This setup was repeated at the Missoula airport (Figure 22). To be given access to the 
airport runway, the team had to submit an experimental plan including risk assessment and 
mitigation around an airport, go through training on airport operations to receive security 
identification display area (SIDA) badges, and supplement the equipment with required airport 
lighting, radios, and signage. 


, a Ms 
Figure 21: Experimental setup at MSU. The camera and asphalt samples were in the temperature-controlled room while data 
collection happened outside the room. Sample and room temperatures were recorded with each image. 
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Figure 22: Experimental setup at Missoula Airport on runway 26. 


Data Analysis Methods 


Data analysis was started by calculating the Road Surface Index value using reflectance 
in the NIR = 1025nm and Red = 389nm where qasphait = 1 with no ice or water present and Qice = - 
1 with no asphalt showing. According to literature, these wavelengths should give best high and 
low reflectance of asphalt and ice. Scripts were written in Python for preliminary image analysis 
using the ice index formula. Each asphalt sample was divided into quadrants and half of those 
quadrants were treated with ice or snow. See Figure 23, quadrants 1,3,5,7 contained 
snow/ice/wet except for September 22" when all quadrants were dry at all temperatures. 


Calibration Panel 


Sample B 


Figure 23: True color image from the Pika L camera of asphalt samples divided into quadrants for analysis. Water can be seen 
in quadrants I and 3 giving off a blueish hue. Colors of the numbers coincide with the graphs under the results section. 


A Python script was written for this task, “Ice_Index_Sample_Mean.py”. This code is 
available in GitHub; https://github.com/jfowler9/MDT Icy Roads.git The script inputs TIFF 
images by wavelength generated from ArcGIS, makes an array of the pixel values, deletes the 
pixel values associated with the calibration panel, then averages ice index values by quadrant of 
each sample. The averages represent both time series and hence room temperature values. 

Another analysis method using a common method for processing and interpreting 
hyperspectral data with principal component analysis (PCA) (Dem$ar et al., 2012) was also 
explored. Due to the size of the datacubes (~1GB each), PCA appeared to be a promising 
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analysis technique. Each layer of the datacube is comprised of pixels of reflectance values for 
each part of the image in one wavelength band. For example, the Resonon Pika L hyperspectral 
camera has 320 wavelength bands resulting in a data cube that is 320 layers deep. Because the 
features in an image may have very similar reflectance values in a number of wavelength bands, 
it is useful to “compress” the information into principal components that explain the majority of 
the variation in the data. These components are linear combinations of the wavelength bands. As 
a simple example, snow may be easily identifiable using 20% of the information from the 400 
nm band, 60% of 420 nm, etc. The Resonon Pika cameras are accompanied by the software 
package Spectronon, one product of which is PCA. ArcGIS also has a PCA product. It was 
decided to build a standalone product that could perform PCA without having to rely on a 
proprietary piece of software and be open-source so a custom Python script was written. This 
code ingests the .bip and .bil files from the hyperspectral cameras directly, converts them to 
raster files, and then to Python arrays. Once in an array format, the machine learning scikit-learn 
libraries are used to perform the PCA. 


SECTION ITI 
Results 
MSU Sub-Zero Laboratory Results — Hyperspectral Imagery 


Results indicate promise for detecting ice versus dry asphalt using UM's Pika L but 
discerning black ice from other ice is not promising using the ice index formula with this current 
dataset as seen in the plots below (Figure 24); 
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Figure 24: plot of ice index value vs time on dry asphalt in all quadrants (variable) for September 22. B; plot of ice index value 
vs time with wet or iced quadrants in 1, 3, 5, and 7 for September 23. Temperature starts at 22 °F (-5.5 °C) - 24 °F (-4.4 °C) each 
day rising to 40 °F (4.4 °C) before descending again from 1400 — 1800. 
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Note immediately that there are no negative ice index values so the basic equation does not hold. 
However, quadrants 2,4, and 6 maintain fairly consistent ice index values from day to day 
staying between .8 and .9. Quadrant 8 is more variable even though it was dry each day. Some of 
this variation calls into question the calibration technique as the only differences between days 
beyond treatment for wet/ice was possible lighting setup. But, the lighting did not change 
throughout a single day hence the oscillation between 14:00 — 18:00 for quadrant 8 is peculiar. 
Ice index values in quadrants 1,3, and 5 decrease from dry to ice/wet as theorized. But, the 
values actually decrease as the water goes through its phase change from ice to liquid 
remembering temperatures are rising from freezing at 7:30 to 40 °F (4.4 °C) by 12:30. Then the 
liquid refreezes between 1400 — 1800. There is also the issue with quadrant 7 with its values 
increasing between dry and ice/wet treatment. Since quadrant 7 and 8 are aligned vertically as 
seen in figure 23 this could be a calibration issue in the software and/or experimental setup. 
Additionally, with such a small range in total index values (.6 - .9) propagating uncertainties 
within the measurements would need to be conducted to confirm the findings are outside the 
margin of error and are statistically significant. 

The Pika L 320 (Figure 25) showed great promise for discerning the difference between 
water and dry asphalt using supervised classification via ArcGIS. Ultimately though, this camera 
requires too much lighting to be of use in operations. 


Figure 25: Supervised classified image from ArcGIS of the Pika L 320 camera depicting water in quadrants I and 3. 


Based on the supervised classification results in figure 25 and, as indicated above, using 
an ice index equation for analysis was not valid, so the processes of exploring the application 
using the python based PCA technique was used with UM’s Pika L. Below are some results from 
the PCA analysis using python (Figures 26 — 28): 
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Figure 26: PCA image for Pika L. Dry asphalt at 21 °F (-6.1 °C). Dark specks are irregularities in the asphalt. 
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Figure 27: PCA image for Pika L at 32 °F (0 °C). Areas of ice may be seen in quadrants 1, 3, 5, and 7 as a darker orange. 
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Figure 28: Pika L at 32 °F (0 °C). 100th wavelength band. Ice visible in quadrants 1, 3, 5, and 7. 


As can be seen in the above figures above this python based PCA technique is capable of 
identifying ice (Figures 27 and 28) in comparison to dry asphalt (Figure 26). Further PCA 
analysis was conducted on the Pike L 320 to double check its possible results using wavelengths 
farther into the NIR (Figures 29 — 36). 
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Figure 29: Pika L 320 camera at 25 °F (-3.8 °C). Dry asphalt. Image in 100th wavelength band. 
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Figure 30: Pika L 320 at 25 °F (-3.9 °C). Dry asphalt PCA. 
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Figure 31: Pika L 320 at 23 °F (-5 °C). Ice shown in 100th wavelength band. 


400 


150 200 250 300 


100 


500 


31 


100 
200 
300 | 
400 


Figure 32: Pika L 320 at 23 °F (-5 °C). Ice visible in PCA. 
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Figure 33: Pika L 320 at 36 °F (2.2 °C). Water visible. 100th wavelength band. 


32 


00 T- 
00s 
si0- 
00r 
os0- 
sz 0- 
OOE 
000 
sz0 = 
oso 
Oot 
slo 
ooT i) 


Figure 34: Pika L 320 at 36 °F (2.2 °C). PCA - water visible in lower left. 
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°C). 100th wavelength band - snow visible. 


Figure 35: Pika L 320 at 26 °F (-3.3 
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Figure 36: Pika L 320 at 23 °F (-5 °C). PCA - snow visible 
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As can be seen in the above figures this python based PCA technique is capable of 
identifying ice (Figures 31 and 32), water (Figures 33 and 34), and snow (Figure 35 and 36) in 
comparison to dry asphalt (Figures 29 and 30). This additional analysis confirmed both UM’s 
Pike L and Resonon’s Pika L 320 are capable of discerning the difference between dry, wet, and 
ice conditions on asphalt, and both could potentially discern the difference between black ice and 
dry asphalt within a controlled laboratory setting using a visual inspection for training a dataset 
to be used in a machine learning model. The Pika L 320 delineates differences better than UM’s 
Pika L but the lighting requirements and weight of that camera disqualify it as useful in UAS 
operations. 


Airport Results — Hyperspectral Imagery 


Figure 37: RGB photo of tread. Area I is ice, Area 2 has an area of black ice, and Area 3 is dry asphalt. 


34 


Based on findings from the MSU sub-zero laboratory analysis, an attempt was made to 
calculate ice index values from the airport setup. Two primary sample areas were established. 
One was called “tread” from vehicle tracks on runway 26. This area had a combination of ice 
(whitish in color), dry asphalt, and black ice (Figure 37). The second area was called “block”. 
This area was shoveled by hand as runway 26 does not receive regular plowing or surface 
treatment while it is closed. This area had a combination of ice (whitish in color), dry asphalt, 
wet, and snow (Figure 39). Below are results from the airport PCA analysis using UM’s Pika L 
camera (Figure 38). 
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Figure 38: Pika L PCA of tread. Snow identifiable in purple as does the calibration puck. Ice/asphalt appears darkest orange. 


As can be seen in Figure 38, the difference between ice, black ice, and dry asphalt are not 
discernable. The same results were found from the analysis of the “block” area. 


Airport Results — Spectrometer 


To investigate band resolution, UM’s Autonomous Aerial Systems Office (AASO) 
deployed an ASD FieldSpec 3 spectroradiometer. The spectrometer has higher band resolution 
(3nm vs 4nm) than the Pika L camera. ASD intensity vs spectral wavelength was graphed (total 
range of 350—2500 nm with a resolution of 1.4 nm for the region 350-1000 nm and 2 nm for the 
region 1001-2500 nm) (Figure 40). Distinguishing features appear between 1000 - 1100nm and 
1800 - 1900nm that better identify spectral differences between snow, ice, wet, and dry 
pavement in comparison to the hyperspectral data. Although, the ASD also does not distinguish 
black ice when it was used on the airport tread data. Also, more work would need to be done to 
create a working model for operational use of the ASD spectrometer including an evaluation of 
its operational feasibility on a UAS considering its weight and size which is in a backpack as 
seen in figure 39 right side. 
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Area 5 is snow. Right: data collection using spectrometer. 


Icy Road_ Project 
ASD_Data_ Block KMSO 
—s — == —< = —= 
Calibration.asd Dry.asd Ice.asd Snow.asd Wet.asd Wet_Dry.asd 


a2 


400 600 800 1000 1200 1400 1600 1800 2000 2200 2400 
Wavelength 


Figure 40: ASD results for the block area at KMSO. 


UM Campus Results — UAS Flight 


To accomplish this last objective, UM’s Pika L camera was mounted onto a Matrice 600 
with a Gremsey H16 gimbal for a flight over a UM parking lot. As noted in the image (Figure 
21, left image), the camera was not able to point at nadir due to a malfunction in the camera 
gimbal. Repair was not possible during the data collection period of this study. But, with this 
image a noticeable decrease in data resolution on the UAS can be seen and a subsequent 
accuracy decrease compared to the airport PCA analysis was found (Figure 41, right image). 
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Figure 41: Left; RGB imagery from Matrice 600 Right; Hyperspectral imagery of Pika L data on Matrice 600 UAS. Note the ice 


cannot be decoupled from the tree shadows 


Possible Future Directions 


PCA is not only a popular and widely used dimensionality reduction technique of 
hyperspectral data, but it also improves classification accuracy (Bajwa et al, 2009). The next 
steps in the project could be to create training data sets (although much more imagery is needed) 
to train a machine learning model to detect ice, snow, wet, and dry pavement automatically. But, 
from this initial analysis, it is not promising that black ice is easily identifiable with UM’s 
hyperspectral camera and the expense associated with collecting enough data for training is not 
cost effective. Additionally, the work to implement a machine learning model that may or may 
not be able to discern black ice is an expensive proposition. Currently this method to confirm 
black ice in a sample is a visual inspection combined with a tap test of the identified black ice 
area and dry asphalt to confirm the black ice acoustically. Once confirmed, the black ice area 
was marked in an image. This process would need to be repeated thousands of times to create a 
training dataset. For this process to have a better chance for success it is suggested to use a 
camera with higher band resolution. 
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